

//GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG
// EdgePrcs.cpp		Implement file
//
// Author:			Cheng-en Guo
//
// Created Time:	07-15-99
// Revised Time:
//
// Descriptions:	This file is to define the class to do edge processing.
//					Including filter, edge detection, build the edge map, and segment by edges
//
// Member variables:
//
// Functions:
//
//GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG

#include "CgEdgePrcs.h"
#include "opencv2\highgui\highgui.hpp"

void CEdgePrcs::Initialize()
{
	m_nColor = 0;// gray level
	char *c = new char[512];
	strcpy( c,
	   "1,2,3,4,5,6,7,8,0,\
		7,8,1,2,3,0,0,0,0,\
		8,1,2,3,4,0,0,0,0,\
		1,2,3,4,5,0,0,0,0,\
		2,3,4,5,6,0,0,0,0,\
		3,4,5,6,7,0,0,0,0,\
		4,5,6,7,8,0,0,0,0,\
		5,6,7,8,1,0,0,0,0,\
		6,7,8,1,2,0,0,0,0," );
	char *temp = c;
	int i;
	for(  i = 0; i < 9; i++ )
	{
		for( int j = 0; j < 9; j++ )
		{
			sscanf( temp, "%d,", &(m_ppnSearchDirection[i][j]) );
			temp += 2;
		}
		temp += 2;
	}
	delete[] c;
	m_pGPointNeighbor[0].x = 1;
	m_pGPointNeighbor[0].y = 0;
	m_pGPointNeighbor[1].x = 1;
	m_pGPointNeighbor[1].y = 1;
	m_pGPointNeighbor[2].x = 0;
	m_pGPointNeighbor[2].y = 1;
	m_pGPointNeighbor[3].x = -1;
	m_pGPointNeighbor[3].y = 1;
	m_pGPointNeighbor[4].x = -1;
	m_pGPointNeighbor[4].y = 0;
	m_pGPointNeighbor[5].x = -1;
	m_pGPointNeighbor[5].y = -1;
	m_pGPointNeighbor[6].x = 0;
	m_pGPointNeighbor[6].y = -1;
	m_pGPointNeighbor[7].x = 1;
	m_pGPointNeighbor[7].y = -1;

	m_nRawImageWidth = 0;
	m_nRawImageHeight = 0;
	m_dSigma = 0.5;
	for( i = 0; i < LEVEL_NUM; i++ )
	{
		m_pdPercentThreshHi[i] = 1 - i/(double)LEVEL_NUM - 0.1;
		m_pdPercentThreshLo[i] = 1 - i/(double)LEVEL_NUM - 0.2;
		m_pdStrengthThreshHi[i] = 0;
		m_pdStrengthThreshLo[i] = 0;
	}
	m_pdPercentThreshHi[0] = 0.60;
	m_pdPercentThreshLo[0] = 0.30;
/*	m_pdPercentThreshHi[1] = 0.00;
	m_pdPercentThreshLo[1] = 0.00;
	m_pdPercentThreshHi[2] = 0.00;
	m_pdPercentThreshLo[2] = 0.00;
	m_pdPercentThreshHi[3] = 0.00;
	m_pdPercentThreshLo[3] = 0.00;
	m_pdPercentThreshHi[4] = 0.00;
	m_pdPercentThreshLo[4] = 0.0;
*/
	m_nMinContourLength = 1;
	m_nMinContourLength = 6;
	m_pContour = NULL;
	m_pCurContour = NULL;

}//void CEdgePrcs::Initialize()

CEdgePrcs::CEdgePrcs()
{
	Initialize();
}

CEdgePrcs::~CEdgePrcs()
{
	if( m_pContour!= NULL )
	{
		delete[] m_pContour;
		m_pContour = NULL;
	}
}

void CEdgePrcs::SetDimension(int nHeight, int nWidth)
{
	m_nRawImageWidth = nWidth;
	m_nRawImageHeight = nHeight;

	InitializeOtherVariables();
}


int CEdgePrcs::InitializeOtherVariables()
{
	m_FMatrixEdgeStrength.SetDimension(m_nRawImageHeight, m_nRawImageWidth);// the float matrix of the edge strength
	m_FMatrixEdgeX.SetDimension(m_nRawImageHeight, m_nRawImageWidth);// the float matrix of the edge strength in x direction
	m_FMatrixEdgeY.SetDimension(m_nRawImageHeight, m_nRawImageWidth);// the float matrix of the edge strength in y direction
	m_FMatrixEdgeAngle.SetDimension(m_nRawImageHeight, m_nRawImageWidth);// the float matrix of the edge angle
	m_PMatrixContour.SetDimension(m_nRawImageHeight, m_nRawImageWidth);
	m_BMatrixLocalMaxFlag.SetDimension(m_nRawImageHeight, m_nRawImageWidth);
	for( int i = 0; i < LEVEL_NUM; i++ )
	{
		m_pBMatrixLocalMax[i].SetDimension(m_nRawImageHeight, m_nRawImageWidth);
		m_pBMatrixLocalMax0[i].SetDimension(m_nRawImageHeight, m_nRawImageWidth);
	}

	return 0;
}

int CEdgePrcs::EdgeDetectionWithEigenvector(Raster_Color<uchar> *praster_color)
{
	double	factor = 180 / PI;


	// Zhuowen Tu's edge detector
	if (praster_color != NULL)
	{
		Raster<double>		Ix, Iy, ras_canny;
		MCMC_EdgeLinking	op_edge_linking;

		if (m_nColor)
		{
			MCMC_EdgeOp_Color edge_op_color;

			edge_op_color = *praster_color;

			// canny edge detector
			ras_canny = edge_op_color.Canny(m_dSigma, 0.12, &Ix, &Iy, true);

		}
		else
		{
			Raster<uchar> ras_intensity;
			MCMC_EdgeOp   edge_op;

			praster_color->GrayImage(ras_intensity);

			edge_op = ras_intensity;

			// canny edge detector
			ras_canny = edge_op.Canny(m_dSigma, 0.1, &Ix, &Iy);
		}
		
		
		op_edge_linking.ObtainAtomicRegions(ras_canny, Ix, Iy, m_dSigma);

		for (int y = 0; y < m_nRawImageHeight; y++)
		{
			for (int x = 0; x < m_nRawImageWidth; x++)
			{
				//m_FMatrixEdgeStrength.Data(y,x) = (float)sqrt(Ix(y,x)*Ix(y,x)+Iy(y,x)*Iy(y,x));
				if (op_edge_linking.ras_edges(y, x) > 0)
					m_FMatrixEdgeStrength(y, x) = (float)max(ras_canny(y, x), 15.0);

				float angle = (float)(atan2(Iy(y, x), Ix(y, x))*factor);

				if (angle >= 0)
					m_FMatrixEdgeAngle(y, x) = angle;
				else
					m_FMatrixEdgeAngle(y, x) = 180 + angle;

				m_FMatrixEdgeX(y, x) = (float)(sqrt(Ix(y, x)*Ix(y, x) + Iy(y, x)*Iy(y, x))*sin(angle / factor));
				m_FMatrixEdgeY(y, x) = (float)(sqrt(Ix(y, x)*Ix(y, x) + Iy(y, x)*Iy(y, x))*cos(angle / factor));
			}
		}

	}
	return 0;
}//int CEdgePrcs::EdgeDetectionWithEigenvector()


int CEdgePrcs::SetPecentThreshHi( int i, double* pd )
{// A bug, no check the decreasing order
	if( i <= 0 || i > LEVEL_NUM )
		return 0;
	for( int j = 0; j < i; j++ )
	{
		m_pdPercentThreshHi[j] = pd[j];
	}
	return i;
}

int CEdgePrcs::SetPecentThreshLo( int i, double* pd )
{// A bug, no check the decreasing order
	if( i <= 0 || i > LEVEL_NUM )
		return 0;
	for( int j = 0; j < i; j++ )
	{
		m_pdPercentThreshLo[j] = pd[j];
	}
	return i;
}
int CEdgePrcs::SetStrengthThreshHi( int i, double* pd )
{// A bug, no check the decreasing order
	if( i <= 0 || i > LEVEL_NUM )
		return 0;
	for( int j = 0; j < i; j++ )
	{
		m_pdStrengthThreshHi[j] = pd[j];
	}
	return i;
}
int CEdgePrcs::SetStrengthThreshLo( int i, double* pd )
{// A bug, no check the decreasing order
	if( i <= 0 || i > LEVEL_NUM )
		return 0;
	for( int j = 0; j < i; j++ )
	{
		m_pdStrengthThreshLo[j] = pd[j];
	}
	return i;
}

double CEdgePrcs::FindPercentThreshold( mcmcMatrix<float> *pFM, double fTh )
{
	if( fTh >= 1.0 )
		return fTh;
	int nWidth = pFM->cols();
	int nHeight = pFM->rows();
#define HISNUM 1000
	double pdHis[HISNUM];
	memset( pdHis, 0, sizeof(double)*HISNUM );
// get min max values
	float fMax = pFM->Data(0, 0);
	float fMin = fMax;
	int y;
	for( y = 0; y < nHeight; y++ )
	{
		for( int x = 0; x < nWidth; x++ )
		{
			if( fMax < pFM->Data(y, x) )
				fMax = pFM->Data(y,x);
			if( fMin > pFM->Data(y,x) )
				fMin = pFM->Data(y,x);
		}
	}
// divide [min, max] into HISNUM equal intervals
	for( y = 0; y < nHeight; y++ )
	{
		for( int x = 0; x < nWidth; x++ )
		{
			int nInterval = (int)( 0.5 + (pFM->Data(y,x) - fMin) * HISNUM / (fMax-fMin) );
			if( nInterval ==  HISNUM )
				nInterval =  HISNUM - 1;
			pdHis[nInterval]++;
		}
	}
	double dMarker = nWidth * nHeight * fTh;
	double dSum = 0;
	int i;
	for( i = 0; i < HISNUM && dSum < dMarker; i++ )
	{
		dSum += pdHis[i];
	}
	double th = (float)(fMin + (i-1) * (fMax-fMin) / HISNUM);
	return th;
}//double CEdgePrcs::FindPercentThreshold( mcmcMatrix<float> *pFM, double fTh )

int CEdgePrcs::GetStrengthThreshFromPercentThresh()
{
	for( int i = 0; i < LEVEL_NUM; i++ )
	{// A bug, if the percent threshs are not set to LEVEL_NUM pairs or not in the decreasing order
		m_pdStrengthThreshHi[i] = FindPercentThreshold( &m_FMatrixEdgeStrength, m_pdPercentThreshHi[i] );
		m_pdStrengthThreshLo[i] = FindPercentThreshold( &m_FMatrixEdgeStrength, m_pdPercentThreshLo[i] );
	}
	return 0;
}//int CEdgePrcs::GetStrengthThreshFromPercentThresh()

int CEdgePrcs::GetEdgeLocalMax()
{
// set local max matrix
	int l = 0;
	m_BMatrixLocalMaxFlag.InitValue( 255 );
	for( l = 0; l < LEVEL_NUM; l++ )
	{
		m_pBMatrixLocalMax[l].InitValue( 0 );
		m_pBMatrixLocalMax0[l].InitValue( 0 );
	}
	l = 0;
	for( int y = 0; y < m_nRawImageHeight; y++ )
	{
		for( int x = 0; x < m_nRawImageWidth; x++ )
		{
			if( EdgeLocalMax( x, y, &m_FMatrixEdgeStrength, &m_FMatrixEdgeX, &m_FMatrixEdgeY, m_pdStrengthThreshLo[l] ) )
			{
//				for( l = 0; l < LEVEL_NUM - 1; l++ )
//				{
//					if( m_FMatrixEdgeStrength.Data(y,x) > m_pdStrengthThreshHi[l] )
//						break;
//				}
//				if( l < LEVEL_NUM )// find
//				{
					m_BMatrixLocalMaxFlag(y,x) = l;
//					for( int ll = l; ll < LEVEL_NUM; ll++ )
//						m_pBMatrixLocalMax[ll].Data(y,x) = 255;
						m_pBMatrixLocalMax[l](y,x) = 255;
//						m_pBMatrixLocalMax0[l].Data(y,x) = 255;
//				}
			}// if local max
		}//for( int x = 0; x < m_nRawImageWidth; x++ )
	}//	for( int y = 0; y < m_nRawImageHeight; y++ )
	return 0;
}//int CEdgePrcs::GetEdgeLocalMax()



BOOL CEdgePrcs::EdgeLocalMax( int ax, int ay, mcmcMatrix<float> *pFMEdgeStrength, mcmcMatrix<float> *pFMEdgeX, mcmcMatrix<float> *pFMEdgeY, double fTh )
{
	int nWidth = pFMEdgeStrength->cols();
	int nHeight = pFMEdgeStrength->rows();
	if( ay <= 0 || ay >= nHeight - 1 || ax <= 0 || ax >= nWidth - 1 )
		return FALSE;
	float mag = pFMEdgeStrength->Data(ay,ax);
	if( mag <= fTh )
		return FALSE;
// edgex edgey give tangent to edge,  use perpenticular
//	float x = pFMEdgeX->Data(ay,ax);
//	float y = -pFMEdgeY->Data(ay,ax);
	float x = pFMEdgeY->Data(ay,ax);
	float y = -pFMEdgeX->Data(ay,ax);

/*  float x = (float)2*( pFMEdgeStrength->ppData[ay][ax+1] -  pFMEdgeStrength->ppData[ay][ax-1] )
				+ ( pFMEdgeStrength->ppData[ay+1][ax+1] -  pFMEdgeStrength->ppData[ay+1][ax-1] )
				+ ( pFMEdgeStrength->ppData[ay-1][ax+1] -  pFMEdgeStrength->ppData[ay-1][ax-1] );
	float y = (float)2*(pFMEdgeStrength->ppData[ay+1][ax] -  pFMEdgeStrength->ppData[ay-1][ax])
				+ (pFMEdgeStrength->ppData[ay+1][ax-1] -  pFMEdgeStrength->ppData[ay-1][ax-1])
				+ (pFMEdgeStrength->ppData[ay+1][ax+1] -  pFMEdgeStrength->ppData[ay-1][ax+1]);
*/
	//	if( x == 0 && y== 0 )
//		x = y = (float)sqrt( 2 );
//		return false;
	float absx = (float)fabs( x );
	float absy = (float)fabs( y );
	float maxc, minc, a, b, c, d;
// four cases
	if( absx < absy )
	{
		maxc = absy;
		minc = absx;
		b = pFMEdgeStrength->Data(ay-1, ax);
		c = pFMEdgeStrength->Data(ay+1, ax);
		// *c*
		// ***
		// *b*
	}
	else //if( absx >= absy )
	{
		maxc = absx;
		minc = absy;
		b = pFMEdgeStrength->Data(ay, ax-1);
		c = pFMEdgeStrength->Data(ay, ax+1);
		// ***
		// b*c
		// ***
	}
/*	else
	{// add by cguo
		maxc = absx;
		minc = absy;
		if( pFMEdgeStrength->ppData[ay-1][ax] + pFMEdgeStrength->ppData[ay+1][ax]
			< pFMEdgeStrength->ppData[ay][ax-1] + pFMEdgeStrength->ppData[ay][ax+1] )
		{
			b = pFMEdgeStrength->ppData[ay-1][ax];
			c = pFMEdgeStrength->ppData[ay+1][ax];
		// *c*
		// ***
		// *b*
		}
		else
		{
			b = pFMEdgeStrength->ppData[ay][ax-1];
			c = pFMEdgeStrength->ppData[ay][ax+1];
		// ***
		// b*c
		// ***
		}
	}
*/
	if( x * y > 0 )// same sign
	{ // case 1 or 3
		a = pFMEdgeStrength->Data(ay-1, ax-1);
		d = pFMEdgeStrength->Data(ay+1, ax+1);
		// **d
		// ***
		// a**
	}
	else
	{
		if( absx < absy )
		{ // case 2 - steep
			a = pFMEdgeStrength->Data(ay-1, ax+1);
			d = pFMEdgeStrength->Data(ay+1, ax-1);
		// d**
		// ***
		// **a
		}
		else
		{ // case 4 shallow
			a = pFMEdgeStrength->Data(ay+1, ax-1);
			d = pFMEdgeStrength->Data(ay-1, ax+1);
		// a**
		// ***
		// **d
		}
	}
/*
	float fDelt;
//	fDelt = (float)-0.0001;
	fDelt = 0;
	float adjusted = (maxc+minc) * mag;
	bool bJudge1 = ( adjusted - ( maxc * b + minc * a ) > fDelt );
	bool bJudge2 = ( adjusted - ( maxc * c + minc * d ) > fDelt );
	bool bJudge = bJudge1 && bJudge2;
	if(bJudge)//(  )
		return TRUE;
	else
		return FALSE;
		*/
	float adjusted = (maxc + minc) * mag;
	if( adjusted > maxc*b + minc*a && adjusted > maxc*c + minc*d )
		return true;
	return false;
}//BOOL CEdgePrcs::EdgeLocalMax( int ax, int ay, mcmcMatrix<float> *pFMEdgeStrength, mcmcMatrix<float> *pFMEdgeX, mcmcMatrix<float> *pFMEdgeY, double fTh )


int CEdgePrcs::GetDirectionFromAngle( float fAngle )
{
	if( fAngle >= 157.5 || fAngle <= -157.5 )
		return 5;
	else if( fAngle >= 112.5 )
		return 4;
	else if( fAngle >= 57.5 )
		return 3;
	else if( fAngle >= 22.5 )
		return 2;
	else if( fAngle >= -22.5 )
		return 1;
	else if( fAngle >= -77.5 )
		return 8;
	else if( fAngle >= -112.5 )
		return 7;
	else
		return 6;
}//int CEdgePrcs::GetDirectionFromAngle( float fAngle )

float CEdgePrcs::GetAngleFromDirection( float fAngle, int nDir )
{
	if( fAngle < 0 )
	{
		::MessageBox( NULL, "Angle is negative", "Anlge error", MB_OK );
		return 0;
	}
	switch ( nDir )
	{
	case 2: case 3: case 4:
		return fAngle;
		break;
	case 6: case 7: case 8:
		return fAngle - 180;
		break;
	case 1:
		if( fAngle <= 90 )
			return fAngle;
		else
			return fAngle - 180;
		break;
	case 5:
		if( fAngle >= 90 )
			return fAngle;
		else
			return fAngle - 180;
		break;
	}
	return 0;
}//float CEdgePrcs::GetAngleFromDirection( float fAngle, int nDir )


BOOL CEdgePrcs::Direction2Vector( int nDir, CFVector *pV )
{
	CFVector v(2);
	if( pV->m_pData == NULL )
	{
		pV->Create( 2 );
	}
	else if( pV->m_nDim != 2 )
	{
		pV->ClearData();
		pV->Create( 2 );
	}
	switch ( nDir )
	{
	case 0:
		pV->m_pData[0] = 0;
		pV->m_pData[1] = 0;
		break;
	case 1:
		pV->m_pData[0] = 1;
		pV->m_pData[1] = 0;
		break;
	case 2:
		pV->m_pData[0] = 1;
		pV->m_pData[1] = 1;
		break;
	case 3:
		pV->m_pData[0] = 0;
		pV->m_pData[1] = 1;
		break;
	case 4:
		pV->m_pData[0] = -1;
		pV->m_pData[1] = 1;
		break;
	case 5:
		pV->m_pData[0] = -1;
		pV->m_pData[1] = 0;
		break;
	case 6:
		pV->m_pData[0] = -1;
		pV->m_pData[1] = -1;
		break;
	case 7:
		pV->m_pData[0] = 0;
		pV->m_pData[1] = -1;
		break;
	case 8:
		pV->m_pData[0] = 1;
		pV->m_pData[1] = -1;
		break;
	default:
		pV->m_pData[0] = 0;
		pV->m_pData[1] = 0;
		break;
	}
	return TRUE;
}//BOOL CEdgePrcs::Direction2Vector( int nDir, CFVector *pV )

CGPoint CEdgePrcs::Direction2Point( int nDir )
{
	CGPoint p;
	switch ( nDir )
	{
	case 0:
		p.x = 0;
		p.y = 0;
		break;
	case 1:
		p.x = 1;
		p.y = 0;
		break;
	case 2:
		p.x = 1;
		p.y = 1;
		break;
	case 3:
		p.x = 0;
		p.y = 1;
		break;
	case 4:
		p.x = -1;
		p.y = 1;
		break;
	case 5:
		p.x = -1;
		p.y = 0;
		break;
	case 6:
		p.x = -1;
		p.y = -1;
		break;
	case 7:
		p.x = 0;
		p.y = -1;
		break;
	case 8:
		p.x = 1;
		p.y = -1;
		break;
	default:
		p.x = 0;
		p.y = 0;
		break;
	}
	return p;
}//CGPoint CEdgePrcs::Direction2Point( int nDir )

float CEdgePrcs::Compatible( CFVector &v, CFVector &v1, CFVector &v2 )
{
#define SPL 0.25881904514 // sin(15 degree): pay 15 degrees per unit distance
	float fPenalty = (float)HUGE_PENALTY;
	float fAngle1 = v1.GetAngle();
	float fAngle2 = v2.GetAngle();
	float fLength  = v.GetLength();
	float fLength1 = v1.GetLength();
	float fLength2 = v2.GetLength();
//	if( fLength < 0.000001 )//|| fLength1 == 0 || fLength2 == 0 )
//		return fPenalty;
	fPenalty = (float)(sin( fabs( fAngle1-fAngle2) * PI / 180 ) + SPL * fLength);
	return fPenalty;
}//float CEdgePrcs::Compatible( CFVector &v, CFVector &v1, CFVector &v2 )

void CEdgePrcs::SetPointer2ContourMatrix( CContour *pContour, SPCMatrix *pPCMatrix, BOOL flag )
{
	GPOSITION pos = pContour->m_PointList.GetHeadPosition();
	CGPoint p;
	while ( pos != NULL )
	{
		p = pContour->m_PointList.GetNext( pos );
		if( flag )
			pPCMatrix->Data(p.y, p.x) = pContour;
		else
			pPCMatrix->Data(p.y, p.x) = NULL;
	}
}//void CEdgePrcs::SetPointer2ContourMatrix( CContour *pContour, SPCMatrix *pPCMatrix, BOOL flag )

float CEdgePrcs::CalculateContourAverageEdgeStrength( CContour *pCon, mcmcMatrix<float> *pFMEdgeSt )
{
	float fStrength = 0;
	GPOSITION pos = pCon->m_PointList.GetHeadPosition();
	CGPoint p;
	while ( pos != NULL )
	{
		p = pCon->m_PointList.GetNext( pos );
		fStrength += pFMEdgeSt->Data(p.y, p.x);
	}
	return (pCon->m_fStrength = fStrength / pCon->m_PointList.GetCount());
}//float CEdgePrcs::CalculateContourAverageEdgeStrength( CContour *pCon, mcmcMatrix<float> *pFMEdgeSt )



void CEdgePrcs::FindNeighborMovingPoints( CGPoint p, int &n, CGPoint &p1, CGPoint &p2, CGPoint &p3 )
{
	p1 = p;
	if( abs( p.x ) != abs( p.y ) )
	{
		n = 1;
		if( abs( p.x ) > abs( p.y ) )
			p1.x += Sign( p1.x );
		else
			p1.y += Sign( p1.y );
	}
	else
	{
		n = 3;
		p2 = p;
		p3 = p;
		p1.x += Sign( p1.x );
		p2.x += Sign( p2.x );
		p2.y += Sign( p2.y );
		p3.y += Sign( p3.y );
	}
}//void CEdgePrcs::FindNeighborMovingPoints( CGPoint p, int &n, CGPoint &p1, CGPoint &p2, CGPoint &p3 )

float CEdgePrcs::CalculateContinuityStrength( CContour *pC1, int &nEnd1, CContour *pC2, int &nEnd2, CGPoint p )
{
	int nLength1 = pC1->m_PointList.GetCount();
	int nLength2 = pC2->m_PointList.GetCount();

	float fFactor[10];
	fFactor[0] = 1;
	fFactor[1] = 0;
	fFactor[2] = 0;
	float fSumFactor = 0;
	int i;
	for( i = 0; i < 3; i++ )
		fSumFactor += fFactor[i];
	CGPoint p1, p2[3];
	float fAngle1, fAngle2[3];
	float fEdgeStrength1 = (float)fabs(pC1->m_fStrength / 128 );
	float fEdgeStrength2 = (float)fabs(pC2->m_fStrength / 128 );
	if( fEdgeStrength1 > 1 )
		fEdgeStrength1 = 1;
	if( fEdgeStrength2 > 1 )
		fEdgeStrength2 = 1;
	if( nEnd1 == 0 )
	{
		p1 = pC1->m_PointList.GetHead();
		fAngle1 = pC1->m_fAngleHead;
	}
	else
	{
		p1 = pC1->m_PointList.GetTail();
		fAngle1 = pC1->m_fAngleTail;
	}
	p2[0] = pC2->m_PointList.GetHead();
	p2[1] = pC2->m_PointList.GetTail();
	p2[2] = p;
	int nPoints4Compare = 3;
	if( p2[2] == p2[0] || p2[2] == p2[1] )
		nPoints4Compare = 2;
	fAngle2[0] = pC2->m_fAngleHead;
	fAngle2[1] = pC2->m_fAngleTail;
	fAngle2[2] = m_FMatrixEdgeAngle.Data(p.y, p.x);
	float fCntStrength[3];
	int nLevel1 = pC1->m_nLevel;
	int nLevel2 = pC2->m_nLevel;
	float fMinStrength = 0;//-(float)(m_nRawImageWidth + m_nRawImageHeight);
	for( i = 0; i < nPoints4Compare; i++ )
	{
		float dis = (float)( abs(p1.x - p2[i].x) + abs(p1.y - p2[i].y) );
		if( dis == 0 )
			dis = fMinStrength;
		else
			dis = (float)(1.0 / dis);
//		float fLinkAngle = atan2( p2[i].x - p1.x, p2[i].y - p1.y )*180/PI;
		float angle = 1//360
			/ (float)(fabs( fabs( fAngle1 - fAngle2[i])-180) + 1);
		float egstrth = ( fEdgeStrength1 + fEdgeStrength2 ) / 2;
		fCntStrength[i] = ( dis*fFactor[0] + angle*fFactor[1] + egstrth*fFactor[2] ) / fSumFactor;
//		fCntStrength[i] = fCntStrength[i] * ( 10 - nLevel1 - nLevel2 ) / 10;
	}
//	fCntStrength[2] /= 2;
	if( pC1 == m_pContourListWholeImage[pC1->m_nLevel].m_pContourHead
		|| pC2 == m_pContourListWholeImage[pC2->m_nLevel].m_pContourHead )
	{
		fCntStrength[2] /= 2;
	}
	else
		fCntStrength[2] = fCntStrength[2]*2/3;
//	for( i = 0; i < nPoints4Compare; i++ )
//	{
//		CContour LinkContour( p1, p2[i] );
//		if( ContourCrossSegLine( LinkContour, m_pBMatrixLocalMax0[nLevel1] ) )
//			fCntStrength[i] = fMinStrength;
//	}
	nEnd2 = 0;
	float fStrength = fCntStrength[0];
	for( i = 1; i < nPoints4Compare; i++ )
	{
		if( fCntStrength[i] >= fCntStrength[nEnd2] )
		{
			nEnd2 = i;
//			fStrength = fCntStrength[nEnd2];
		}
	}
	return fCntStrength[nEnd2];
}//float CEdgePrcs::CalculateContinuityStrength( CContour *pC1, int &nEnd1, CContour *pC2, int &nEnd2, CGPoint p )



bool CEdgePrcs::PointsNeighborJudge( CGPoint p1, CGPoint p2 )
{
	if( abs( p1.x - p2.x ) <= 1 && abs( p1.y - p2.y ) <= 1 )
		return true;
	else
		return false;
}//bool CEdgePrcs::PointsNeighborJudge( CGPoint p1, CGPoint p2 )


int CEdgePrcs::SetBMatrixFlagOfPointList( mcmcMatrix<BYTE>& BMatrixFlag, CGPointList& PointList, BYTE bValue )
{
	GPOSITION pos = PointList.GetHeadPosition();
	CGPoint p;
	while ( pos != NULL )
	{
		p = PointList.GetNext( pos );
		if( p.x >= 0 && p.x < BMatrixFlag.cols() && p.y >= 0 && p.y < BMatrixFlag.rows() )
			BMatrixFlag(p.y,p.x) = bValue;
	}
	return 0;
}//int CEdgePrcs::SetBMatrixFlagOfPointList( mcmcMatrix<BYTE>& BMatrixFlag, CGPointList& PointList, BYTE bValue )

bool CEdgePrcs::ContourCrossSegLine( CContour& LinkContour, mcmcMatrix<BYTE>& BMatrixOnSegLineFlag )
{
	int nWidth = BMatrixOnSegLineFlag.cols();
	int nHeight = BMatrixOnSegLineFlag.rows();
	GPOSITION pos = LinkContour.m_PointList.GetHeadPosition();
	CGPoint p, p1, p2;
	p = LinkContour.m_PointList.GetNext( pos );// skip the first point
//	if( pos != NULL )
//		p = LinkContour.m_PointList.GetNext( pos );// skip the first point
	while ( pos != NULL )
	{
		if( BMatrixOnSegLineFlag.Data(p.y,p.x) == 255 )
			return true;
		//below deal with diagonal cross
		p1.x = p.x - 1;
		p1.y = p.y;
		p2.x = p.x;
		p2.y = p.y - 1;
		p1 = MovePointIntoRegion( p1, nWidth, nHeight );
		p2 = MovePointIntoRegion( p2, nWidth, nHeight );
		if( BMatrixOnSegLineFlag.Data(p1.y,p1.x) == 255
			&& BMatrixOnSegLineFlag.Data(p2.y,p2.x) == 255 )
			return true;
		p2.y = p.y + 1;
		p2 = MovePointIntoRegion( p2, nWidth, nHeight );
		if( BMatrixOnSegLineFlag.Data(p1.y,p1.x) == 255
			&& BMatrixOnSegLineFlag.Data(p2.y,p2.x) == 255 )
			return true;
		p1.x = p.x + 1;
		p1 = MovePointIntoRegion( p1, nWidth, nHeight );
		if( BMatrixOnSegLineFlag.Data(p1.y,p1.x) == 255
			&& BMatrixOnSegLineFlag.Data(p2.y,p2.x) == 255 )
			return true;
		p2.y = p.y - 1;
		p2 = MovePointIntoRegion( p2, nWidth, nHeight );
		if( BMatrixOnSegLineFlag.Data(p1.y,p1.x) == 255
			&& BMatrixOnSegLineFlag.Data(p2.y,p2.x) == 255 )
			return true;
		p = LinkContour.m_PointList.GetNext( pos );
	}
	return false;
}//bool CEdgePrcs::ContourCrossSegLine( CContour& LinkContour, mcmcMatrix<BYTE>& BMatrixOnSegLineFlag )


int CEdgePrcs::SegmentRegions(Raster_Color<uchar> &raster_color)
{
	int nReturn = 0;

	nReturn = EdgeDetectionWithEigenvector(&raster_color);

	nReturn = GetStrengthThreshFromPercentThresh();
	if( nReturn != 0 )
		return nReturn;
	//set over segment boundary matrice
	for( int i = 0; i < LEVEL_NUM; i++ )
	{
		m_pBMatrixOverSegmentLines[i].SetDimension(m_nRawImageHeight, m_nRawImageWidth);
		m_pBMatrixOverSegmentLines[i].InitValue( 0 );
	}
	for( int l = 0; l < 1; l++ )
	{
		nReturn = BuildContourList4OverSegment( l );
		if( nReturn != 0 )
			return nReturn;
///			nReturn = BuildContourMap4OverSegment( l, 1 );//no neighbor for no extend contour

		nReturn = De_freeTheEndsOfAllContours( l );

		if( nReturn != 0 )
			return nReturn;
	}
	if( nReturn != 0 )
		return nReturn;


	m_MatrixRegionLabels.SetDimension(m_nRawImageHeight, m_nRawImageWidth);
	m_MatrixRegionLabels.InitValue(-1);

	LabelRegions( m_pBMatrixLocalMax0[0], m_MatrixRegionLabels, m_nRegionNum );

	return nReturn;

}//int CEdgePrcs::OverSegmentRegions()

int CEdgePrcs::De_freeTheEndsOfAllContours( int nLevel )
{
	int nReturn = 0;
	CContour* pCurContour = m_pContourListWholeImage[nLevel].m_pContourHead;
	while( pCurContour != NULL )
	{//for each contour in this level, de-free its' ends
	// extension this contour
	// what is considered is in the whole image, not in the specialized region
		for( int i = 0; i < 2; i++ )
		{//consider two direction, 0-head, 1-tail
			nReturn = De_freeOneEndOfOneContour( pCurContour, NULL, i, nLevel );
		}//for( int dir = 0; dir < 2; dir++ )
		pCurContour = pCurContour->m_pNext;
	}//while( pCurContour != NULL )
	return nReturn;
}//int CEdgePrcs::De_freeTheEndsOfAllContoursInEdgeMap()

int CEdgePrcs::De_freeOneEndOfOneContour( CContour* pCurContour, CContour* pContourAlreadyLinked, int i, int nLevel )
{
	int nReturn = 0;
	int nWidth = m_nRawImageWidth;
	int nHeight = m_nRawImageHeight;

	CContourNeighbor CntNbr;
	CContour* pNbrContour;
	int nEndNbrContour;
	CGPoint p1;//the end of pCurContour
	CGPoint p2;//the end of pNbrContour
			while( pCurContour->m_bExtendHeadTail[i] )
			{//still needs extension
				if( ( (pCurContour->m_PointList.m_pNodeHead->m_data.y) == 108
//					&& (pCurContour->m_PointList.m_pNodeHead->m_data.y) == 188
						)
					||
					( (pCurContour->m_PointList.m_pNodeTail->m_data.y) == 108
//						&& (pCurContour->m_PointList.m_pNodeTail->m_data.y) == 188
						) )
				{
					int test = 0;
				}
				GetNeighborsOfAContourAtOneEnd(pCurContour, nLevel, 1, i );//1-means if meet boudary then stop
				pCurContour->SortNeighbors( i );
				int nNbrNum = pCurContour->GetNeighborNum( i );
				if( nNbrNum <= 0 )
				{// no neighbor
//					DirectExtendContourUntilNoneFree( pCurContour, nLevel, i );
//					pCurContour->m_bExtendHeadTail[i] = false;
//					pCurContour->m_bExtendHeadTail[1-i] = false;
#ifdef _DEBUG
//					::MessageBox( NULL, "The neighbor number of one contour is 0", "Error", MB_OK );
#endif//_DEBUG
					SetBMatrixFlagOfPointList( m_pBMatrixLocalMax0[nLevel], pCurContour->m_PointList, 0 );
					break;
				}
				BOOL bExtend = false;
				float fExtendNbrStrength;
				for( int j = 0; j < nNbrNum; j++ )
				{
					CntNbr = pCurContour->GetNeighbor( i, j );
					nEndNbrContour = CntNbr.m_nEnd;
					pNbrContour = CntNbr.m_pContour;
					if( pNbrContour == pContourAlreadyLinked )
					{
							if( !bExtend )
							{
								bExtend = true;
								fExtendNbrStrength = CntNbr.m_fContinuityStrength;
								pCurContour->m_bExtendHeadTail[i] = false;
							}
							continue;
					}
					if( bExtend )
					{// extend other neighbors
//						break;
						if( CntNbr.m_fContinuityStrength <= fExtendNbrStrength * 3/4 )
							break;
						if( nEndNbrContour == 2 )
							continue;// the link point on the neighbor is inside point, not end
						p1 = pCurContour->GetEndPoint( i );
						p2 = pNbrContour->GetEndPoint( nEndNbrContour );
						if( PointsNeighborJudge( p1, p2 ) )
						{
//							pNbrContour->m_bExtendHeadTail[nEndNbrContour] = false;
							continue;
						}
						CContour LinkContour( p1, p2 );
						if( ContourCrossSegLine( LinkContour, m_pBMatrixLocalMax0[nLevel] ) )
						{//link contour cross other contour
							continue;
						}
						else
						{//it is allowed to link
//							CombineContour( pNbrContour, nEndNbrContour, &LinkContour, 0, nLevel, 0 );//don't delete link contour
							SetBMatrixFlagOfPointList( m_pBMatrixLocalMax0[nLevel], LinkContour.m_PointList, 255 );
//							pNbrContour->m_bExtendHeadTail[nEndNbrContour] = false;
						}
						continue;
					}
					if( pNbrContour == pCurContour )
					{// next neighbor is the original contour,
						p1 = pCurContour->GetEndPoint( i );
						if( nEndNbrContour != 2 )// end point
							p2 = pNbrContour->GetEndPoint( nEndNbrContour );
						else
							p2 = CntNbr.m_Point;// inside point
						if( !PointsNeighborJudge( p1, p2 ) )
						{// p1 and p2 are not neighbor
							CContour LinkContour( p1, p2 );
							if( ContourCrossSegLine( LinkContour, m_pBMatrixLocalMax0[nLevel] ) )
							{//link contour cross other contour
								continue;
							}
							else
							{//it is allowed to link
//								CombineContour( pCurContour, i, &LinkContour, 0, nLevel, 0 );//don't delete link contour
								SetBMatrixFlagOfPointList( m_pBMatrixLocalMax0[nLevel], LinkContour.m_PointList, 255 );
//								if( nEndNbrContour == 2 )// cross the inside of itself
//									CutContour( pCurContour, p2, 0, nLevel );
							}
						}
						if( i + nEndNbrContour == 1 )
						{//the line is closed
							pCurContour->m_bExtendHeadTail[0] = false;
							pCurContour->m_bExtendHeadTail[1] = false;
						}
						else if( nEndNbrContour == 2 )
						{// reach the inside point of the contour
							pCurContour->m_bExtendHeadTail[i] = false;
						}
						else
						{// reach the same end of itself
							pCurContour->m_bExtendHeadTail[i] = false;
#ifdef _DEBUG
//							::MessageBox( NULL, "One contour reaches itself at the same end", "Error", MB_OK );
#endif
						}
						bExtend = true;
						fExtendNbrStrength = CntNbr.m_fContinuityStrength;
						continue;
					}//if( pNbrContour == pCurContour )
					// else, the neighbor is not itself
					p1 = pCurContour->GetEndPoint( i );
					if( nEndNbrContour != 2 )//
						p2 = pNbrContour->GetEndPoint( nEndNbrContour );
					else
						p2 = CntNbr.m_Point;
					if( !PointsNeighborJudge( p1, p2 ) )
					{
						CContour LinkContour( p1, p2 );
						if( ContourCrossSegLine( LinkContour, m_pBMatrixLocalMax0[nLevel] ) )
						{//link contour cross other contour
							if( nEndNbrContour == 2 )
								continue;// next contour
							p2 = CntNbr.m_Point;
							CContour LinkContour2( p1, p2 );
							if( ContourCrossSegLine( LinkContour2, m_pBMatrixLocalMax0[nLevel] ) )
								continue;
							LinkContour = LinkContour2;
						}
//						else
						{//it is allowed to link
//							CombineContour( pCurContour, i, &LinkContour, 0, nLevel, 0 );//don't delete link contour
							SetBMatrixFlagOfPointList( m_pBMatrixLocalMax0[nLevel], LinkContour.m_PointList, 255 );
							pCurContour->m_bExtendHeadTail[i] = false;
							if( nEndNbrContour != 2 )
								De_freeOneEndOfOneContour( pNbrContour, pCurContour, nEndNbrContour, nLevel );
//								CombineContour( pCurContour, i, pNbrContour, , nLevel, 1 );
	//							pNbrContour->m_bExtendHeadTail[nEndNbrContour] = false;
//							else
//								CutContour( pNbrContour, p2, 0, nLevel );
							bExtend = true;
							fExtendNbrStrength = CntNbr.m_fContinuityStrength;
							continue;//next neighbor contour
						}
					}
					else
					{//p1, p2 neighbor
							if( nEndNbrContour != 2 )
								De_freeOneEndOfOneContour( pNbrContour, pCurContour, nEndNbrContour, nLevel );
//								CombineContour( pCurContour, i, pNbrContour, nEndNbrContour, nLevel, 1 );
	//							pNbrContour->m_bExtendHeadTail[nEndNbrContour] = false;

	//						else
	//						{
//								CutContour( pNbrContour, p2, 0, nLevel );
//								pCurContour->m_bExtendHeadTail[i] = false;
	//						}
							pCurContour->m_bExtendHeadTail[i] = false;
							bExtend = true;
							fExtendNbrStrength = CntNbr.m_fContinuityStrength;
							continue;//next neighbor contour
					}
				}//for( int j = 0; j < nNbrNum; j++ )
				if( !bExtend )
				{
#ifdef _DEBUG
//					::MessageBox( NULL, "One contour is not de-freed", "Error", MB_OK );
#endif
//					DirectExtendContourUntilNoneFree( pCurContour, nLevel, i );
					pCurContour->m_bExtendHeadTail[i] = false;
				}
			}//while( pCurContour->m_bExtendHeadTail[i] )
			return 0;
}//int CEdgePrcs::De_freeOneEndOfOneContour( CContour* pCurContour, CContour* pContourAlreadyLinked, int nEnd, int nLevel )

int CEdgePrcs::BuildContourList4OverSegment( int nLevel )
{
	int nReturn = 0;
// Get edge local max first
	GetEdgeLocalMax();
// below get the contours
	m_PMatrixContour.InitValue( NULL );// this matrix saves the pointers to the contour which the current point belongs to
	CContour *pNewContour = NULL;
	CGPoint p;
	int l = 0;
// below set the boader as a contour
	pNewContour = new CContour;
	int x, y;
	for( x = 0; x < m_nRawImageWidth; x++ )
	{
		p.x = x;
		p.y = 0;
		pNewContour->m_PointList.AddTail( p );
		m_BMatrixLocalMaxFlag(p.y,p.x) = nLevel;
		m_PMatrixContour(p.y,p.x) = pNewContour;
		m_FMatrixEdgeAngle(p.y,p.x) = 0;
	}
	for( y = 1; y < m_nRawImageHeight; y++ )
	{
		p.x = m_nRawImageWidth - 1;
		p.y = y;
		pNewContour->m_PointList.AddTail( p );
		pNewContour->m_PointList.AddTail( p );
		m_BMatrixLocalMaxFlag(p.y,p.x) = nLevel;
		m_PMatrixContour(p.y,p.x) = pNewContour;
		m_FMatrixEdgeAngle(p.y,p.x) = 90;
	}
	for( x = m_nRawImageWidth - 2; x >= 0 ; x-- )
	{
		p.x = x;
		p.y = m_nRawImageHeight - 1;
		pNewContour->m_PointList.AddTail( p );
		pNewContour->m_PointList.AddTail( p );
		m_BMatrixLocalMaxFlag(p.y,p.x) = nLevel;
		m_PMatrixContour(p.y,p.x) = pNewContour;
		m_FMatrixEdgeAngle(p.y,p.x) = 180;
	}
	for( y = m_nRawImageHeight - 2; y > 0; y-- )
	{
		p.x = 0;
		p.y = y;
		pNewContour->m_PointList.AddTail( p );
		pNewContour->m_PointList.AddTail( p );
		m_BMatrixLocalMaxFlag(p.y,p.x) = nLevel;
		m_PMatrixContour(p.y,p.x) = pNewContour;
		m_FMatrixEdgeAngle(p.y,p.x) = 90;
	}
	pNewContour->m_bExtendHeadTail[0] = false;
	pNewContour->m_bExtendHeadTail[1] = false;
	pNewContour->m_nLevel = nLevel;
	m_pContourListWholeImage[nLevel].AddTail( pNewContour );
	SetBMatrixFlagOfPointList( m_pBMatrixLocalMax0[nLevel], pNewContour->m_PointList, 255 );
	for( y = 0; y < m_nRawImageHeight; y++ )
	{
		for( x = 0; x < m_nRawImageWidth; x++ )
		{//Search every point in the image,
			l = m_BMatrixLocalMaxFlag.Data(y,x);
			if( m_PMatrixContour.Data(y,x) != NULL ||  l != nLevel )
				continue;
			if( nLevel >= 1 && m_pBMatrixOverSegmentLines[nLevel].Data(y,x) == 255 )
				continue;
			else
			{
			 //if it is local max and does not belong to other contour
			// and at the same level as nLevel, then find a new contour
				pNewContour = new CContour;
				p.x = x;
				p.y = y;
				GetContourAtOnePoint4OverSegment(
					nLevel,
					p,
					pNewContour,
					m_pdStrengthThreshHi[l],
					m_pdStrengthThreshLo[l],
					0// nFlag - 0: trace only if finding the point at the same level
						//      1: trace continuously unless meet the bounarary
					);
				pNewContour->m_nLevel = l;
				if( m_nMinContourLength > 1 )
				{
					if( pNewContour->m_PointList.GetCount() >= m_nMinContourLength )
					{
//						SetPointer2ContourMatrix( pNewContour, &m_PMatrixContour, TRUE );
//						CalculateContourAverageEdgeStrength( pNewContour, &m_FMatrixEdgeStrength );
						m_pContourListWholeImage[l].AddTail( pNewContour );
						SetBMatrixFlagOfPointList( m_pBMatrixLocalMax0[nLevel], pNewContour->m_PointList, 255 );
					}
					else
					{
						SetPointer2ContourMatrix( pNewContour, &m_PMatrixContour, FALSE );
//						m_pCurContour--;
						delete pNewContour;
					}
				}
				else
				{
//						CalculateContourAverageEdgeStrength( pNewContour, &m_FMatrixEdgeStrength );
						m_pContourListWholeImage[l].AddTail( pNewContour );
						SetBMatrixFlagOfPointList( m_pBMatrixLocalMax0[nLevel], pNewContour->m_PointList, 255 );
				}
			}//if( m_PMatrixContour.Data(y,x) == NULL && EdgeLocalMax(...)
		}//for( int x = 0; x < m_nRawImageWidth; x++ )
	}//for( int y = 0; y < m_nRawImageHeight; y++ )
	return nReturn;
}//int CEdgePrcs::BuildContourList4OverSegment( int l )

void CEdgePrcs::GetContourAtOnePoint4OverSegment(
					int nLevel,
					CGPoint p,
					CContour *pNewContour,
					double fThreshHi,
					double fThreshLo,
					int nFlag
					)
{
	int x = p.x;
	int y = p.y;
	pNewContour->m_PointList.AddHead( p );
	m_PMatrixContour.Data(y,x) = pNewContour;
	pNewContour->m_fAngleHead = m_FMatrixEdgeAngle.Data(y,x);
	pNewContour->m_fAngleTail = pNewContour->m_fAngleHead - 180;
	if( pNewContour->m_fAngleTail == 0 )
		pNewContour->m_fAngleTail = -180;
	for( int i = 0; i <2; i++ )
	{
		int nFlag1 = 0;
		pNewContour->m_nLastMoveDirection = 0;
		while ( nFlag1 == 0 && pNewContour->m_PointList.GetCount() <= m_nMaxContourLength
			&& pNewContour->m_bExtendHeadTail[i] )
		{
			nFlag1 = TraceContour4OverSegment( nLevel, i, pNewContour,fThreshHi, fThreshLo, nFlag );
		}
	}
}//void CEdgePrcs::GetContourAtOnePoint4OverSegment(...)

int CEdgePrcs::TraceContour4OverSegment( int nLevel, int end, CContour *pContour, double fThreshHi, double fThreshLo, int nFlag )
{
	int nReturn = 0;//0-continue tracing, 1-stop tracing
	int nWidth = m_FMatrixEdgeStrength.cols();
	int nHeight = m_FMatrixEdgeStrength.rows();
	CGPoint p0, p1;//p0 is the point which the search is from, p1 is candidate extension point
	if( pContour->m_nLastMoveDirection == 0 )
	{// no last move direction, get this direction from the tail point angle
		if( end == 1 )
			pContour->m_nLastMoveDirection = GetDirectionFromAngle( pContour->m_fAngleTail );
		else // extend head
			pContour->m_nLastMoveDirection = GetDirectionFromAngle( pContour->m_fAngleHead );
	}
	if( end == 1 )
		p0 = pContour->m_PointList.GetTail();
	else
		p0 = pContour->m_PointList.GetHead();
	if( ( p0.x == 48 && p0.y == 200 )
		|| ( p0.x == 49 && p0.y == 200 )
		|| ( p0.x == 51 && p0.y == 200 )
		)
		int test = 0;
	int d = pContour->m_nLastMoveDirection;// get move direction
	p1.x = p0.x;
	p1.y = p0.y;
	if( p0.x < 1 || p0.x > nWidth - 2 || p0.y < 1 || p0.y > nHeight -  2 )
	{// reach the boundary, mark it
		pContour->m_nHeadAndTailProperty[end] = 2;
		pContour->m_bExtendHeadTail[end] = false;
		nReturn = 1;//stop tracing
		return nReturn;
	}
	int nBestDir = -1;
	float fBestCompat = (float)HUGE_PENALTY;
	float fCompat = fBestCompat;
	CGPoint dP;
	CFVector v;
	CFVector v1(2);
	CFVector v2(2);
	int sd;// search point direction
	for( int i = 0; i < 9; i++ )
	{
		sd = m_ppnSearchDirection[d][i];// search point direction
		if( sd == 0 )
			break;
		dP = Direction2Point( sd );
		p1.x = p0.x + dP.x;
		p1.y = p0.y + dP.y;
		if( p1.x < 1 || p1.x > nWidth - 2 || p1.y < 1 || p1.y > nHeight -  2 )
		{// reach the boundary, mark it
			pContour->m_nHeadAndTailProperty[end] = 2;
			nReturn = 1;//stop tracing
			nBestDir = sd;// search point direction
			pContour->m_bExtendHeadTail[end] = false;
			break;
		}
		if( m_BMatrixLocalMaxFlag.Data(p1.y,p1.x) != nLevel )
			continue;//for searching only in given level
		if( m_PMatrixContour.Data(p1.y,p1.x) != NULL )
		{//meet other contour, stop, divide the met contour
			CContour* pMetContour = m_PMatrixContour.Data(p1.y,p1.x);
			if( pMetContour == pContour )
			{// meet itself
				int nCloseLength = pContour->GetLengthFromInsidePoint2OneEnd( end, p1 );
				if( nCloseLength >= 10*m_dSigma )
				{//large close area
					if( p1 != pContour->m_PointList.GetHeadOrTail( end ) )
					{// meet the inside of itself, cut it, and keep the cut part to a new contour
						CutContour( pContour, p1, end, nLevel );
					}
					//else no need cut
					pContour->m_nHeadAndTailProperty[end] = 0;//reach itself
					pContour->m_nHeadAndTailProperty[1-end] = 0;
	//				pContour->m_bExtendHeadTail[end] = false;
	///				pContour->m_bExtendHeadTail[1-end] = false;
					nReturn = 1;
					nBestDir = sd;// search point direction
					return nReturn;
				}
				else
					continue;
			}
			else
			{// meet other contour
				if( p1 == pMetContour->m_PointList.GetHead() )
				{//meet the head of that contour
					CombineContour( pContour, end, pMetContour, 0, nLevel, 1 );
					return 1;
				}
				else if( p1 == pMetContour->m_PointList.GetTail() )
				{//meet the tail of that contour
					CombineContour( pContour, end, pMetContour, 1, nLevel, 1 );
					return 1;
				}
				else
				{//meet inside of that contour
					CutContour( pMetContour, p1, end, nLevel );
					pContour->m_nHeadAndTailProperty[end] = 1;
					nReturn = 1;
					nBestDir = sd;// search point direction
					return nReturn;
				}
			}//else meet other contour
		}//if( m_PMatrixContour.Data(p1.y,p1.x) != NULL )
		Direction2Vector( sd, &v );
		v1.m_pData[0] = m_FMatrixEdgeX.Data(p0.y,p0.x);
		v1.m_pData[1] = m_FMatrixEdgeY.Data(p0.y,p0.x);
		v2.m_pData[0] = m_FMatrixEdgeX.Data(p1.y,p1.x);
		v2.m_pData[1] = m_FMatrixEdgeY.Data(p1.y,p1.x);
		fCompat = Compatible( v, v1, v2 );
		if( nBestDir == -1 || fBestCompat > fCompat )
		{// only searching in the same level
			nBestDir = sd;
			fBestCompat = fCompat;
		}
	}
	// nFlag - 0: trace only if finding the point at the same level
	//      1: tracing could be at different level
	//      2: tracing in direction as the last even though local maximum is not there,
	if( nBestDir == -1 )// not find the next search point in the same level
	{//stop tracing
			nReturn = 1;
			return nReturn;
	}// else, find tracing point at the same level, the the direction is "nBestDir"

	pContour->m_nLastMoveDirection = nBestDir;
	dP = Direction2Point( nBestDir );
	p1.x = p0.x + dP.x;
	p1.y = p0.y + dP.y;

	// new valid point, add it to the contour
	if( end == 1 )
	{
		pContour->m_PointList.AddTail( p1 );
		pContour->m_fAngleTail = GetAngleFromDirection( m_FMatrixEdgeAngle.Data(p1.y,p1.x), nBestDir );
	}
	else
	{
		pContour->m_PointList.AddHead( p1 );
		pContour->m_fAngleHead = GetAngleFromDirection( m_FMatrixEdgeAngle.Data(p1.y,p1.x), nBestDir );
	}
	m_PMatrixContour.Data(p1.y,p1.x) = pContour;
	return nReturn;
}//int CEdgePrcs::TraceContour4OverSegment( int nLevel, int end, CContour *pContour, double fThreshHi, double fThreshLo, int nFlag )

int CEdgePrcs::CutContour( CContour* pContour, CGPoint p1, int end, int nLevel )
{

	int nReturn = 0;
	CContour* pCutContour = new CContour;
	CGPoint p;
	if( end == 0 )//head
	{//keep the head for pContour, cut the tail for pCutcontour, search from tail
		do
		{
			p = pContour->m_PointList.RemoveTail();
			pCutContour->m_PointList.AddHead( p );
			m_PMatrixContour.Data(p.y,p.x) = pCutContour;
		}while( p != p1 );
		pCutContour->m_bExtendHeadTail[0] = false;
		pCutContour->m_bExtendHeadTail[1] = pContour->m_bExtendHeadTail[1];
		pCutContour->m_nLevel = nLevel;
		pCutContour->m_nHeadAndTailProperty[0] = 1;//reach other contour
		pCutContour->m_nHeadAndTailProperty[1] = pContour->m_nHeadAndTailProperty[1];
		memcpy( pCutContour->m_TailPtNbrArray, pContour->m_TailPtNbrArray, sizeof(CContourNeighbor)*pContour->m_nTailNbrNum );
		pContour->m_nTailNbrNum = 0;
	}
	else
	{//keep the tail for pContour, cut the head for pCutcontour, search from head
		do
		{
			p = pContour->m_PointList.RemoveHead();
			pCutContour->m_PointList.AddTail( p );
			m_PMatrixContour.Data(p.y,p.x) = pCutContour;
		}while( p != p1 );
		pCutContour->m_bExtendHeadTail[0] = pContour->m_bExtendHeadTail[0];
		pCutContour->m_bExtendHeadTail[1] = false;
		pCutContour->m_nLevel = nLevel;
		pCutContour->m_nHeadAndTailProperty[1] = 1;//reach other contour
		pCutContour->m_nHeadAndTailProperty[0] = pContour->m_nHeadAndTailProperty[0];
		memcpy( pCutContour->m_HeadPtNbrArray, pContour->m_HeadPtNbrArray, sizeof(CContourNeighbor)*pContour->m_nHeadNbrNum );
		pContour->m_nHeadNbrNum = 0;
	}
// add the cut contour to the whole image
	m_pContourListWholeImage[nLevel].AddTail( pCutContour );

	return nReturn;
}//int CEdgePrcs::CutContour( CContour* pContour, CGPoint p1, int end, int nLevel )

int CEdgePrcs::CombineContour( CContour* pC1, int end1, CContour* pC2, int end2, int nLevel, int nDelete )
{
	int nReturn = 0;

	GPOSITION pos;
	CGPoint p;
	if( end1 == 0 )//head
	{
		if( end2 == 0 )
		{//combine pC2'head to pC1's head
			pos = pC2->m_PointList.GetHeadPosition();
			while( pos != NULL )
			{
				p = pC2->m_PointList.GetNext( pos );
				pC1->m_PointList.AddHead( p );
				m_PMatrixContour.Data(p.y,p.x) = pC1;
			}
			memcpy( pC1->m_HeadPtNbrArray, pC2->m_TailPtNbrArray, sizeof(CContourNeighbor)*pC2->m_nTailNbrNum );
			pC1->m_nHeadNbrNum = pC2->m_nTailNbrNum;
			pC1->m_bExtendHeadTail[0] = pC2->m_bExtendHeadTail[1];
		}
		else
		{//combine pC2's tail to pC1's head
			pos = pC2->m_PointList.GetTailPosition();
			while( pos != NULL )
			{
				p = pC2->m_PointList.GetPrev( pos );
				pC1->m_PointList.AddHead( p );
				m_PMatrixContour.Data(p.y,p.x) = pC1;
			}
			memcpy( pC1->m_HeadPtNbrArray, pC2->m_HeadPtNbrArray, sizeof(CContourNeighbor)*pC2->m_nHeadNbrNum );
			pC1->m_nHeadNbrNum = pC2->m_nHeadNbrNum;
			pC1->m_bExtendHeadTail[0] = pC2->m_bExtendHeadTail[0];
		}
	}
	else
	{//pC1's tail
		if( end2 == 0 )
		{//combine pC2'head to pC1's tail
			pos = pC2->m_PointList.GetHeadPosition();
			while( pos != NULL )
			{
				p = pC2->m_PointList.GetNext( pos );
				pC1->m_PointList.AddTail( p );
				m_PMatrixContour.Data(p.y,p.x) = pC1;
			}
			memcpy( pC1->m_TailPtNbrArray, pC2->m_TailPtNbrArray, sizeof(CContourNeighbor)*pC2->m_nTailNbrNum );
			pC1->m_nTailNbrNum = pC2->m_nTailNbrNum;
			pC1->m_bExtendHeadTail[1] = pC2->m_bExtendHeadTail[1];
		}
		else
		{//combine pC2's tail to pC1's tail
			pos = pC2->m_PointList.GetTailPosition();
			while( pos != NULL )
			{
				p = pC2->m_PointList.GetPrev( pos );
				pC1->m_PointList.AddTail( p );
				m_PMatrixContour.Data(p.y,p.x) = pC1;
			}
			memcpy( pC1->m_TailPtNbrArray, pC2->m_HeadPtNbrArray, sizeof(CContourNeighbor)*pC2->m_nHeadNbrNum );
			pC1->m_nTailNbrNum = pC2->m_nHeadNbrNum;
			pC1->m_bExtendHeadTail[1] = pC2->m_bExtendHeadTail[0];
		}
	}
	if( nDelete == 1 )
		m_pContourListWholeImage[nLevel].Remove( pC2 );

	return nReturn;
}//int CEdgePrcs::CombineContour( CContour* pC1, int end1, CContour* pC2, int end2, int nLevel, int nDelete )

int	CEdgePrcs::BuildContourMap4OverSegment(int nLevel, int nFlag)
// nFlag 0-all possible neighbors, 1-no neighbors if reach boundary or other contour
{
	int nWidth = m_nRawImageWidth;
	int nHeight = m_nRawImageHeight;
// blow get the contour map, which is to find the neighbors of each contour
	CGPoint *pPtArray1 = new CGPoint[(m_nRawImageWidth+m_nRawImageHeight)*2];
	CGPoint *pPtArray2 = new CGPoint[(m_nRawImageWidth+m_nRawImageHeight)*2];
	int nPtArraySize1 = 0;
	int nPtArraySize2 = 0;
	CGPoint p[2];//p[0]-the head point of a contour, p[1]-the tail point of a contour
	CGPoint pt, dP, ptN, p1, p2, p3;
//	for( int l = 0; l < LEVEL_NUM; l++ )
	int l = nLevel;
	{
		CContour *pCurContour = m_pContourListWholeImage[l].m_pContourHead;//get the head contour
		while( pCurContour != NULL )
		{
			if( pCurContour->m_nTheEndMeetOtherContour == 2 && nFlag == 1 )
			{
				pCurContour = pCurContour->m_pNext;// next contour
				continue;//both ends meet other contour
			}
			if( pCurContour->m_nTheEndReachBoundary == 2 && nFlag == 1 )
			{
				pCurContour = pCurContour->m_pNext;// next contour
				continue;//both ends reach boundary
			}
			if( pCurContour->m_bExtendHeadTail[0] == false && pCurContour->m_bExtendHeadTail[1] == false )
			{
				pCurContour = pCurContour->m_pNext;// next contour
				continue;//both ends reach boundary
			}
			p[0] = pCurContour->m_PointList.GetHead();
			p[1] = pCurContour->m_PointList.GetTail();
			for( int i = 0; i < 2; i++ )
			{// look for the neighbors both for the head(i=0) and for the tail(i=1)
				if( pCurContour->m_nTheEndMeetOtherContour == i && nFlag == 1 )
				{//the current end meet other contour
					continue;
				}
				if( pCurContour->m_nTheEndReachBoundary == i && nFlag == 1 )
				{//the current end reach boundary
					continue;
				}
				if( pCurContour->m_bExtendHeadTail[i] == false )
				{
					continue;//both ends reach boundary
				}
				pt = p[i];
				nPtArraySize1 = 8;
				nPtArraySize2 = 0;
				memcpy( pPtArray1, m_pGPointNeighbor, 8*sizeof(CGPoint) );
				// PtArray1 is for the currently being checked points, it is initialized to 8 connection points
				// PtArray2 is for further checked points
				BOOL bContinue = FALSE;
				int nNeighborNum = pCurContour->GetNeighborNum( i );
				CContourNeighbor *pNbrArray;
				if( i == 0 )
				{
					pNbrArray = pCurContour->m_HeadPtNbrArray;
				}
				else
				{
					pNbrArray = pCurContour->m_TailPtNbrArray;
				}
				while ( (nPtArraySize1 > 0 || nPtArraySize2 > 0) && nNeighborNum < MAX_NEIGHBOR_NUM )
				{
					if( nPtArraySize1 <= 0 )
					{
						memcpy( pPtArray1, pPtArray2, sizeof(CGPoint) * nPtArraySize2 );
						nPtArraySize1 = nPtArraySize2;
						nPtArraySize2 = 0;
					}
					dP = pPtArray1[nPtArraySize1-1];
					nPtArraySize1--;
					ptN = pt + dP;
					CContour *pCheckingContour = m_PMatrixContour.Data(ptN.y,ptN.x);
					if( ptN.x <= 0 || ptN.x >= nWidth - 1 || ptN.y <= 0 || ptN.y >= nHeight - 1 )
					{//reach the boudary
						pCheckingContour = &m_ContourRepresentBoudary;//point to the the reservered m_pContour[0]
						break;//won't search anymore
					}
					else
					{// not reach boundary
						if( pCheckingContour == NULL )
						{// no contour at ptN, then add 1 or 3 points to PtArray2
							int n;
							FindNeighborMovingPoints( dP, n, p1, p2, p3 );
							if( n == 1 )
							{
								pPtArray2[nPtArraySize2] = p1;
								nPtArraySize2++;
							}
							else if( n == 3 )
							{
								pPtArray2[nPtArraySize2] = p1;
								nPtArraySize2++;
								pPtArray2[nPtArraySize2] = p2;
								nPtArraySize2++;
								pPtArray2[nPtArraySize2] = p3;
								nPtArraySize2++;
							}
							continue;
						}
						// there is a contour at PtN,
						if( pCheckingContour == pCurContour )
							continue;
						// check if current checking contour is already in the neighbor list
						bContinue = FALSE;
						for( int j = 0; j < nNeighborNum; j++ )
						{
							if( pCheckingContour == pNbrArray[j].m_pContour )
							{//alread in the neighbor list, next point
								bContinue = TRUE;
								break;
							}
						}
						if( bContinue )
							continue;//next point
					}//else not reach boundary
					// Add the checking contour to the neighbor
					CContourNeighbor ConNeighbor;
					ConNeighbor.m_pContour = pCheckingContour;
					int nEndN = 0;
					float fStrength = CalculateContinuityStrength( pCurContour, i, pCheckingContour, nEndN, ptN );
					if( fStrength <= 0 )// too weak to link as a neighbor
						continue;
					ConNeighbor.m_fContinuityStrength = fStrength;
					ConNeighbor.m_nEnd = nEndN;
					pCurContour->AddNeighbor( i, ConNeighbor );
			// also add neighbor for pCheckingContour //no, it is not correct when the neighbor number is limited
//					ConNeighbor.m_pContour = pCurContour;
//					ConNeighbor.m_nEnd = i;
//					pCheckingContour->AddNeighbor( nEndN, ConNeighbor );
				}//while ( PtArray1.GetSize() > 0 || PtArray2.GetSize() > 0 )
			}//for( int i = 0; i < 2; i++ )
			pCurContour = pCurContour->m_pNext;// next contour
		}//while( pCurContour != NULL )
	}
	delete[] pPtArray1;
	delete[] pPtArray2;
	return 0;
}//int	CEdgePrcs::BuildContourMap4OverSegment( int nLevel, int nFlag )

int	CEdgePrcs::GetNeighborsOfAContourAtOneEnd( CContour* pCurContour, int nLevel, int nFlag, int nEnd )
// nFlag 0-all possible neighbors, 1-no neighbors if reach boundary or other contour
{
	int nWidth = m_nRawImageWidth;
	int nHeight = m_nRawImageHeight;
// blow get the contour map, which is to find the neighbors of each contour
	CGPoint *pPtArray1 = new CGPoint[(m_nRawImageWidth+m_nRawImageHeight)*2];
	CGPoint *pPtArray2 = new CGPoint[(m_nRawImageWidth+m_nRawImageHeight)*2];
	int nPtArraySize1 = 0;
	int nPtArraySize2 = 0;
	CGPoint p;//the end point of a contour
	CGPoint pt, dP, ptN, p1, p2, p3;
	int l = nLevel;
//	if( pCurContour->m_PointList.GetCount() <= 0 )
//		return 1;
//	if( pCurContour->m_bExtendHeadTail[nEnd] == false )
//		return 0;
	p = pCurContour->m_PointList.GetHeadOrTail( nEnd);
	pt = p;
	nPtArraySize1 = 8;
	nPtArraySize2 = 0;
	memcpy( pPtArray1, m_pGPointNeighbor, 8*sizeof(CGPoint) );
	// PtArray1 is for the currently being checked points, it is initialized to 8 connection points
	// PtArray2 is for further checked points
	BOOL bContinue = FALSE;
	int nNeighborNum;
//	pCurContour->ClearNeighbor( nEnd );
	CContourNeighbor *pNbrArray;
	if( nEnd == 0 )
		pNbrArray = pCurContour->m_HeadPtNbrArray;
	else
		pNbrArray = pCurContour->m_TailPtNbrArray;
	CContour *pCheckingContour;
	BOOL bBreakWhile = false;
	CContourNeighbor ConNeighbor;
	int nContourLength = pCurContour->m_PointList.GetCount();
	int nSearchRadius = 0;
	int nSearchLimit = (int)(nContourLength + m_dSigma);
	while ( (nPtArraySize1 > 0 || nPtArraySize2 > 0) && pCurContour->GetNeighborNum(nEnd) < MAX_NEIGHBOR_NUM )
	{
		if( pCurContour->m_bExtendHeadTail[1-nEnd] && nSearchRadius > nSearchLimit )
			break;// the other end is free and the search area is too large
		if( nPtArraySize1 <= 0 )
		{
			memcpy( pPtArray1, pPtArray2, sizeof(CGPoint) * nPtArraySize2 );
			nPtArraySize1 = nPtArraySize2;
			nPtArraySize2 = 0;
			nSearchRadius++;
		}
		dP = pPtArray1[nPtArraySize1-1];
		nPtArraySize1--;
		ptN = pt + dP;
		if( ptN.x == 39 && ptN.y == 214 )
			int test = 0;
		if( ptN.x < 0 || ptN.x > nWidth - 1 || ptN.y < 0 || ptN.y > nHeight - 1 )
			continue;
//		{//reach the boudary
/*			pCheckingContour = new CContour;
			ptN = MovePointIntoRegion( ptN, nWidth, nHeight );
			pCheckingContour->m_PointList.AddTail( ptN );
			pCheckingContour->m_nHeadAndTailProperty[0] = 3;
			pCheckingContour->m_nHeadAndTailProperty[1] = 3;
			pCheckingContour->m_nLevel = nLevel;
			pCheckingContour->m_bExtendHeadTail[0] = false;
			pCheckingContour->m_bExtendHeadTail[1] = false;
			m_pContourListWholeImage[nLevel].AddTail( pCheckingContour );
			m_PMatrixContour.ppData[ptN.y][ptN.x] = pCheckingContour;
//			m_pBMatrixLocalMax[nLevel].ppData[ptN.y][ptN.x] = 255;
			m_BMatrixLocalMaxFlag.ppData[ptN.y][ptN.x] = nLevel;
//			bBreakWhile = true;
*/
//		}
		else
		{// not reach boundary
			pCheckingContour = m_PMatrixContour.Data(ptN.y,ptN.x);
			if( pCheckingContour == NULL )
			{// no contour at ptN, then add 1 or 3 points to PtArray2
				int n;
				FindNeighborMovingPoints( dP, n, p1, p2, p3 );
				if( n == 1 )
				{
					pPtArray2[nPtArraySize2] = p1;
					nPtArraySize2++;
				}
				else if( n == 3 )
				{
					pPtArray2[nPtArraySize2] = p1;
					nPtArraySize2++;
					pPtArray2[nPtArraySize2] = p2;
					nPtArraySize2++;
					pPtArray2[nPtArraySize2] = p3;
					nPtArraySize2++;
				}
				continue;//next point
			}
			// there is a contour at PtN,
			if( pCheckingContour == pCurContour )//the same as pCurContour, next point
			{
				if( ptN == pCurContour->m_PointList.GetHeadOrTail( 1 - nEnd ) )
				{
					int nLength = pCurContour->GetLengthFromInsidePoint2OneEnd( nEnd, ptN );
					if( nLength <= 20*m_dSigma || nLength <= 2 )
					{
						int test = 0;
						continue;
					}
				}
				else
					continue;
			}
			// check if current checking contour is already in the neighbor list
			bContinue = FALSE;
			nNeighborNum = pCurContour->GetNeighborNum( nEnd );
			for( int j = 0; j < nNeighborNum; j++ )
			{
				if( pCheckingContour == pNbrArray[j].m_pContour )
				{//alread in the neighbor list, next point
					bContinue = TRUE;
					break;
				}
			}
			if( bContinue )
				continue;//next point
		}//else not reach boundary
		// Add the checking contour to the neighbor
		ConNeighbor.m_pContour = pCheckingContour;
		int nEndN = 0;
		float fStrength = CalculateContinuityStrength( pCurContour, nEnd, pCheckingContour, nEndN, ptN );
//		if( fStrength <= 0 )// too weak to link as a neighbor
//			continue;
		ConNeighbor.m_fContinuityStrength = fStrength;
		ConNeighbor.m_nEnd = nEndN;
//		if( ConNeighbor.m_nEnd == 2 )
		ConNeighbor.m_Point = ptN;
		pCurContour->AddNeighbor( nEnd, ConNeighbor );
		if( bBreakWhile )
			break;
	}//while ( PtArray1.GetSize() > 0 || PtArray2.GetSize() > 0 )
	delete[] pPtArray1;
	delete[] pPtArray2;
	return 0;
}//int	CEdgePrcs::GetNeighborsOfAContourAtOneEnd( CContour* pCurContour, int nLevel, int nFlag, int nEnd )

CGPoint CEdgePrcs::MovePointIntoRegion( CGPoint ptN, int nWidth, int nHeight )
{
	CGPoint pReturn = ptN;
	if( ptN.x < 0 )
		pReturn.x = 0;
	if( ptN.x >= nWidth )
		pReturn.x = nWidth - 1;
	if( ptN.y < 0 )
		pReturn.y = 0;
	if( ptN.y >= nHeight )
		pReturn.y = nHeight - 1;
	return pReturn;
}//CGPoint CEdgePrcs::MovePointIntoRegion( CGPoint ptN, int nWidth, int nHeight )

int CEdgePrcs::DirectExtendContourUntilNoneFree( CContour* pCurContour, int nLevel, int end )
{
	CGPoint p0 = pCurContour->GetEndPoint( end );
	float fAngle;
	if( end == 0 )
		fAngle = pCurContour->m_fAngleHead;
	else
		fAngle = pCurContour->m_fAngleTail;
	return 0;
}

int CEdgePrcs::LabelRegions( mcmcMatrix<BYTE>& BMatrixEdge, mcmcMatrix<int> & MatrixRegionLabels, int& nLabelNo )
{
	nLabelNo = 0;
	CGPoint* pPoints = new CGPoint[m_nRawImageHeight*m_nRawImageWidth];
	for( int y = 0; y < m_nRawImageHeight; y++ )
	{
		for( int x = 0; x < m_nRawImageWidth; x++ )
		{
			if( BMatrixEdge.Data(y,x) != 255
				&& MatrixRegionLabels(y, x) == -1 )
			{
				LabelRegionFromOneSeed( x, y, &BMatrixEdge, &MatrixRegionLabels, nLabelNo, pPoints );
				nLabelNo++;
			}
		}
	}
	delete[] pPoints;
	LabelBoundaries( &m_MatrixRegionLabels);
	return 0;
}


int	CEdgePrcs::LabelRegionFromOneSeed( int x, int y, mcmcMatrix<BYTE>* pBMatrixEdge, mcmcMatrix<int>* pMatrixRegionLabels, int nLabelNo, CGPoint* pPoints )
{
	CGPoint pNbrPoints[4];
	pNbrPoints[0].x = 1;
	pNbrPoints[0].y = 0;
	pNbrPoints[1].x = 0;
	pNbrPoints[1].y = 1;
	pNbrPoints[2].x = -1;
	pNbrPoints[2].y = 0;
	pNbrPoints[3].x = 0;
	pNbrPoints[3].y = -1;

	int nPointsNumMinus1 = 0;
	pPoints[0].x = x;
	pPoints[0].y = y;
	int x1, y1, x2, y2;
	(*pMatrixRegionLabels)(pPoints[nPointsNumMinus1].y, pPoints[nPointsNumMinus1].x) = nLabelNo;
	while( nPointsNumMinus1 >= 0 )
	{
		x1 = pPoints[nPointsNumMinus1].x;
		y1 = pPoints[nPointsNumMinus1].y;
		for( int i = 0; i <4; i++ )
		{
			x2 = x1 + pNbrPoints[i].x;
			y2 = y1 + pNbrPoints[i].y;
			if( x2 >= 0 && x2 < m_nRawImageWidth && y2 >= 0 && y2 < m_nRawImageHeight
				&& pMatrixRegionLabels->Data(y2, x2) == -1
				&& pBMatrixEdge->Data(y2,x2) != 255 )
			{
				(*pMatrixRegionLabels)(y2, x2) = nLabelNo;
				pPoints[nPointsNumMinus1].y = y2;
				pPoints[nPointsNumMinus1].x = x2;
				nPointsNumMinus1++;
			}
		}
		nPointsNumMinus1--;
	}
	return 0;
}//int	CEdgePrcs::LabelRegionFromOneSeed( int x, int y, mcmcMatrix<BYTE>& BMatrixEdge, CG_Matrix & MatrixRegionLabels, int nLabelNo )


void CEdgePrcs::LabelBoundaries( mcmcMatrix<int> *pMatrixRegionLabels)
{
	int nWidth = pMatrixRegionLabels->cols();
	int nHeight = pMatrixRegionLabels->rows();

	mcmcMatrix<int> tempMCMCMatrix;

	tempMCMCMatrix = *pMatrixRegionLabels;

	int nMaxWH = max( nWidth, nHeight );
	int x1, y1;
	bool bContinue = true;//
	while ( bContinue )
	{
		bContinue = false;
		for( int y = 0; y < nHeight; y++ )
		{
			for( int x = 0; x < nWidth; x++ )
			{
				if( tempMCMCMatrix(y,x) < 0 )
				{
					x1 = x + 1;
					y1 = y;//left neighbor
					if( x1 >= 0 && x1 < nWidth && y1 >= 0 && y1 < nHeight
						&& tempMCMCMatrix(y1, x1) >= 0 )
					{
						(*pMatrixRegionLabels)(y, x) = tempMCMCMatrix(y1, x1);
						continue;
					}
					x1 = x;
					y1 = y + 1; //above neighbor
					if( x1 >= 0 && x1 < nWidth && y1 >= 0 && y1 < nHeight
						&& tempMCMCMatrix(y1, x1) >= 0 )
					{
						(*pMatrixRegionLabels)(y, x) = tempMCMCMatrix(y1, x1);
						continue;
					}
					x1 = x - 1;
					y1 = y; // right neighbor
					if( x1 >= 0 && x1 < nWidth && y1 >= 0 && y1 < nHeight
					&& tempMCMCMatrix(y1, x1) >= 0 )
					{
						(*pMatrixRegionLabels)(y, x) = tempMCMCMatrix(y1, x1);
						continue;
					}
					x1 = x;
					y1 = y - 1;
					if( x1 >= 0 && x1 < nWidth && y1 >= 0 && y1 < nHeight
						&& tempMCMCMatrix(y1, x1) >= 0 )
					{
						(*pMatrixRegionLabels)(y, x) = tempMCMCMatrix(y1, x1);
						continue;
					}
					bContinue = true;
				}//if( tempMCMCMatrix.data[y][x] < 0 )
			}//for( int x = 0; x < nWidth; x++ )
		}//for( int y = 0; y < nHeight; y++ )
		
		tempMCMCMatrix = *pMatrixRegionLabels;

	}//while ( bContinue )
}//void CEdgePrcs::LabelBoundaries( CG_Matrix *pMatrixRegionLabels);




int CEdgePrcs::LabelRegions(double dSigma, mcmcMatrix<int> &mx_region_labels, bool bColor, Raster_Color<uchar> *praster_color)
{
	int nMinContourLength = 2;;
	int nMaxContourLength = 100000;
	
	SetDimension(praster_color->rows(), praster_color->cols());

	SetColor(bColor);

	SetSigma( dSigma );
	SetMinContourLength( nMinContourLength );
	SetMaxContourLength( nMaxContourLength );
	double dHighThresh = 0.90;
	double dLowThresh = 0.2;

	SetPecentThreshHi( 1, &dHighThresh );
	SetPecentThreshLo( 1, &dLowThresh );

	SegmentRegions(*praster_color);

	mx_region_labels = m_MatrixRegionLabels;

	return GetRegionNum();
}


